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Beyond our Solar System, aurorae have been inferred from radio observations of 


isolated brown dwarfs’”. Within our Solar System, giant planets have auroral 
emission with signatures across the electromagnetic spectrum including infrared 
emission of H,” and methane. Isolated brown dwarfs with auroral signatures in 

the radio have been searched for corresponding infrared features, but only null 
detections have been reported’. CWISEP J193518.59-154.620.3. (W1935 for short) is 
an isolated brown dwarf with a temperature of approximately 482 K. Here we report 
James Webb Space Telescope observations of strong methane emission from W1935 
at 3.326 um. Atmospheric modelling leads us to conclude that a temperature 
inversion of approximately 300 K centred at 1-10 mbar replicates the feature. This 
represents an atmospheric temperature inversion for a Jupiter-like atmosphere 
without irradiation from a host star. A plausible explanation for the strong inversion 
is heating by auroral processes, although other internal and external dynamical 
processes cannot be ruled out. The best-fitting model rules out the contribution 

of H,* emission, which is prominent in Solar System gas giants. However, this is 
consistent with rapid destruction of H,” at the higher pressure where the W1935 


emission originates’. 


Brown dwarfs area class of object that links planetary and stellar astro- 
physics. They have temperatures between approximately 3,000 and 
250 K and spectral classifications of L, T and Y (refs. 5,6). The Y dwarfs 
are arecent addition to our assortment of known compact objects. 
They are probably the coldest sources that formed through the star for- 
mation process’. These cold objects are directly comparable to Jupiter, 
with the coldest known Y dwarf—WISEJ085510.83-071442.5—having a 
temperature of approximately 250 K (ref. 8), just 100 K warmer than 
Jupiter’. Y dwarfs present an extraordinary observational challenge 
for ground-based telescopes given their intrinsic faintness and need 
for infrared instrumentation”. The James Webb Space Telescope 
(JWST), a space-based 6.5 m infrared observatory, is perfectly suited 
for revolutionizing our understanding of brown dwarfs and in turn 
exo-Jupiter atmospheres”. In this work we report observations of 
two brown dwarfs obtained with JWST cycle 1 Guest Observer (GO) 
programme 2124. We have obtained NIRSpec G395H spectra and 
mid-infrared FIOOOW, F1280W and F1800W photometry from the 
Mid-Infrared Instrument (MIRI) on JWST for the Y dwarfs CWISEP 
J193518.59-154620.3 (W1935 for short) and WISE J222055.31-362817.4 
(W2220 for short). 


We combined all literature data on these two objects alongside the 
JWST data for each source and created absolute spectral energy distri- 
butions (SEDs), which we could compare and contrast. By integrating 
over the SEDs using the open-source code SEDkit (ref. 14), we find that 
the luminosity ratio to the Sun (L,,,/L.) for W2220 and W1935 are identi- 
cal within uncertainties, with values of log(L,,,,/L,) equal to -6.4 + 0.1 
and -6.3 + 0.1, respectively. Neither source has any age indications, 
so we assumed a conservative age range of 4.5 + 4.0 Gyr, which we 
used to semi-empirically calculate the following values for W1935 and 
W2220, respectively: radius in Jupiter radii (R,,,) of 0.95 + 0.14 Rjup and 
0.94 + 0.14 Rjup; effective temperature (Tep) of 482 + 38 Kand480 + 41K; 
surface gravity (log g) of 4.7 + 0.5 dex (both); and mass in Jupiter masses 
(Mjup) Of 6-35 Mj, (both). Given that these objects are indistinguish- 
able in their fundamental parameters, they are excellent objects for 
intercomparisons. 

Spectroscopically, we found that W1935 and W2220 are near clones 
of each other, with both showing clear and strong signatures of CH,, 
CO, CO,, H,O and NH; (Fig. 1). There is one visually striking difference 
between the spectral characteristics of the two sources. Although 
W2220 shows the expected CH, q-branch absorption centred at 
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Fig. 1|JWST NIRSpec G395H spectra for the Y dwarfs W1935 and W2220. 
a, NIRSpec G395H portion of the SED for W1935 (black) compared to that of 


Wavelength (um) 


CH,v; band. b, The residuals between the two spectra. 


W2220 (green). Shading for both sources represents the lo uncertainty onthe 


3.326 um, W1935 shows emission over that same wavelength range 
(see inset for the 3.25-3.45 um area in Fig. 1). 

To model the 3.326 um emission feature in W1935 as well as com- 
pare and contrast with W2220, we used the Brewster retrieval frame- 
work’>"®, which has successfully retrieved the properties of numerous 


Fig. 2|JWST G395H spectrum for W1935 overlaid with the best-fitting models 
with and without temperature inversion. Overlaid in purple and dark cyan 
are the median retrieval models for the source with and without temperature 
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length range covered by our data: H,O, CH,, CO, 


flux. Major opacity sources are labelled. The inset is a zoom-in of the 3.326 um 


brown dwarfs” ”. For our baseline model for both objects, we assumed 
a cloud-free atmosphere. Alongside continuum opacity sources 
described in detail in refs. 15,16, we included the following gases as 


an impact in the wave- 
CO,,NH;,,H,S and PH}. 


inversion, respectively. Data uncertainties are shaded in grey. The 67% confidence 
intervals inthe model posteriors are shaded in dark cyan and purple but are of 
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Fig. 3 | Retrieved thermal profiles for W1935 with and without thermal 
inversion. Median posterior profiles for the models are shown in purple and 
darkcyan, respectively. The 67% confidence intervals for the model posteriors 
areindicated with shading for each. Also plotted with dashed lines are the 
condensation curves (assuming solar composition gas) for possible cloud 
species. 


The power of a retrieval is in its ability to parameterize the temperature- 
pressure (7/P) profile, thus giving an insight into the energy distribu- 
tion within a given atmosphere. Brewster can do this ina flexible way, 
without prescribing that the atmosphere be in radiative-convective 
equilibrium or havea particular slope. This gives the 7/P profile freedom 
to adopt whatever shape is justified by the data, including an inversion 
where the temperature increases with altitude. 

The results of the retrieval verified that the two sources were near 
clones in all abundances (Extended Data Table 2). However, the 7/P 
profiles for the two sources show striking differences. Although the 
best-fitting retrieved 7/P profile for W2220 is consistent with the tem- 
perature decreasing with increasing altitude throughout the atmos- 
phere (Extended Data Fig. 3), as expected for an isolated source, the 
best-fitting retrieval for W1935 displays a temperature inversion of 
approximately 300 K centred at approximately 1-10 mbar (Fig. 3). As 
aresult, alongside our baseline model, we tested a model that forbade 
an inversion. For W2220, the ‘no inversion’ result was indistinguishable 
from the base result. However, for W1935, the ‘no inversion’ retrieval 
could not reproduce the CH, emission feature (Fig. 2). Hence, we con- 
clude that the observed CH, emission arises as a result of thermal inver- 
sion in the stratosphere of this cool, isolated brown dwarf. 

Temperature inversions have been inferred before in substellar 
atmospheres, for both brown dwarfs” and giant exoplanets”, not to 
mention being nearly ubiquitous within the Solar System. The com- 
mon feature of all of these cases is the presence of an irradiating star. 
However, the stratospheres of the gas planets in the Solar System dis- 
play temperatures even higher than can be attributed to irradiation 
alone”, Our result represents a spectacular extension of this gas 
giant phenomenon without any stellar irradiation. 

Much work has been dedicated to understanding the Solar System 
cases of enhanced stratospheric heating. Both external heating by auro- 
ral processes and internal energy transport by vertically propagating 
waves from deeper in the atmosphere are possible mechanisms” ”. The 
latter is a plausible explanation for the thermal inversion modelled for 
W1935. However, one would expect this process to be ubiquitous across 
arange of atmospheres. Given that this is the only non-irradiated exam- 
ple to date, sucha universal mechanism is less likely to be responsible. 


External heating by auroral processes may be a more likely mecha- 
nism. Recent observations by ref. 28 indicate that the bulk of the heat- 
ing inJupiter’s upper atmosphere is driven by redistribution from hot 
auroral polar regions. In addition, alongside methane fluorescence 
from solar pumping, some Jovian methane emission has been tenta- 
tively attributed to heating by auroral processes”. 

Ultracool dwarfs, which comprise the lowest-mass stars and 
substellar-mass objects, have long been surmised to host aurora akin 
to those found inthe giant planets of our Solar System, like Jupiter and 
Saturn. Studies have shown that approximately 5% of ultracool dwarfs 
demonstrate highly circularly polarized, rotationally modulated, radio 
pulses attributed to the electron cyclotron maser instability, which 
is the mechanism responsible for auroral radius emission!*°. Such a 
low-detection rate suggests that any stratospheric heating arising 
from those auroral processes should be similarly rare. 

Asa further test for W1935, we implemented a retrieval that included 
H,",acommon emitter produced by aurorae in giant planets. Interest- 
ingly, the addition did not improve the fit for W1935 and yielded a null 
detection for H,” emission. Although the thermal inversion in W1935 
has a similar overlying column mass to the equivalent region in the 
Jovian atmosphere, the higher surface gravity results in a gas density 
that is approximately 100 times higher. At these densities, the lifetimes 
of H,* ions are much shorter than the typical emission timescale’, so 
their absence is not surprising. 

The detection of CH, in the emission from an object with a mass 
range of 6-35 Mju and a temperature of 482 + 38 K (approximately 
300 K warmer than Jupiter) is enticing. Moreover, the appearance of 
atemperature inversion in an object that lacks an irradiating star com- 
pounds the interest. For Solar System giants with equivalent spectral 
emission and upper atmospheric heating, acontributing factor outside 
of solar pumping is auroral processes linked to nearby moons (lo for 
Jupiter and Enceladus for Saturn). Regardless of what is causing the 
thermal inversion on and consequent methane emission from W1935, 
this source represents an outstanding laboratory for investigating 
linked phenomena that are prominent in our own Solar System. 
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Methods 


Sample 

JWST cycle 1 GO programme 2124 (PI Faherty) is obtaining G395H NIR- 
Spec spectra for 12 cold brown dwarfs that span late-type T and Y dwarf 
classes. Sources were selected by their position in the mid-infrared 
(Spitzer) colour-magnitude diagram (M,, 5 versus [3.6]-[4.5] where 
[3.6] and [4.5] are two Spitzer bands). Four colour bins were roughly 
defined to represent changing temperatures (from 800 to less than 
350 K). Three sources were chosen in each colour bin such that one 
source was close to the median property line for the population’. The 
other two were, respectively, the brightest and faintest in M,, 5, for 
that colour bin. Two of the sources chosen were CWISEP J193518.59- 
154620.3 (W1935; ref. 31) and WISE J222055.31-362817.4 (W2220; 
ref. 32), which are the subject of this paper. At the time of the project 
design, we speculated that the [3.6]—-[4.5] colour would define the 
temperature bin and we did not expect W1935 and W2220 to be so 
comparable. The results of this paper are strong evidence that M,,5)is 
a better temperature binning parameter. Extended Data Fig. 1 shows 
the Spitzer colour-magnitude diagram for cold brown dwarfs and 
highlights the sample for programme 2124 with the positions of W1935 
and W2220 emphasized. 


CWISEP J193518.59-154620.3. CWISEPJ193518.59-154620.3 (W1935 
for short) was first reported in ref. 31 after a concerted effort was 
applied to find cold compact sources within newly analysed Wide- 
field Infrared Survey Explorer (WISE) catalogue data****. The object’s 
discovery was the result of a collaboration between the CatWISE 
team and the citizen science project Backyard Worlds: Planet 9 
(ref. 35). A Spitzer follow-up resulted in a Spitzer [3.6]-[4.5] colour 
of 3.24 + 0.31 mag, making it one of the reddest mid-infrared sources 
known to date. Additional follow-up observations were done to 
obtain Spitzer data with a higher signal-to-noise ratio (S/N) in both 
channels, and the colour was refined to a Spitzer [3.6]-[4.5] colour of 
2.984 + 0.034 mag in ref. 36. Reference 7 reported that the parallax for 
this source was 69.3 + 3.8 mas and that the total proper motion was 
293.4 + 16.3 mas yr“. As noted by ref. 36, this source has a particularly 
low tangential velocity (vtan) compared to other Y dwarfs analysed. 
The estimated spectral type from its photometry and parallax was > Y1, 
and the temperature was predicted to be 367 + 79 K (ref. 7). Reference 37 
obtained Gemini NIRI imaging for the source and reported Mauna Kea 
Observatories (MKO) /= 23.93 + 0.33 mag. 


WISE J222055.31-362817.4. WISE J222055.31-362817.4 (W2220 for 
short) was first reported in ref. 32 after a search of the WISE cata- 
logue’! for cold compact objects. Reference 32 followed this up with 
observations by the Spitzer Space Telescope, Keck/NIRSpec-N3 and 
Keck/NIRSpec-N5S as well as AAT/IRIS2 and SOAR/OSIRIS. The source 
was confirmed as a cold brown dwarf with WISE W1, W2, Spitzer [3.6], 
[4.5], and MKO JH photometry as well as near-infrared spectra. 
References 39,40 obtained a Hubble Space Telescope (HST) grism 
spectrum for the source, which confirmed the YO classification with 
higher S/N data. Follow-up photometric and astrometric imaging was 
done using both space-based (for example, Spitzer) and ground-based 
(for example, the FourStar Infrared Camera installed on Magellan) 
instruments and the most recent and highest S/N trigonometric par- 
allax is 95.5 + 2.1 mas (see refs. 7,40—43 for the astrometric history of 
the source). Reference 7 estimated a temperature of 452 + 88 K for 
this source based on its position on the colour-magnitude diagram. 


The data 

JWST programme 2124 obtained both NIRSpec G395H spectra and 
MIRI FIOOOW, F1280W, and F1800W photometry to fill out the peak 
of the SED and the Rayleigh—Jeans tail of the SED for 12 brown dwarfs. 
NIRSpec data were obtained using the F290LP filter, the G395H grating, 


the S200A1 aperture and the SUB2048 subarray. The resultant wave- 
length coverage ranged from 2.87 to 5.14 um witha resolving power of 
approximately 2,700. Acquisition images were first obtained for each 
target using the WATA method, the CLEAR filter and the NRSRAPID 
readout pattern. W2220 was observed with NIRSpec on 4 November 
2022 with 28 groups per integration, ten integrations per exposure 
and three dithers for a summation of 30 total integrations in 1,356 s 
of exposure time. The recorded time including the overhead for the 
W2220 NIRSpec observation was 1.42 h. W1935 was observed with NIR- 
Spec on17 October 2022 with 46 groups per integration, llintegrations 
per exposure and three dithers for asummation of 33 total integrations 
in 2,417 s of exposure time. The recorded time including the overhead 
for the W1935 NIRSpec observation was 1.76 h. 

MIRI photometry was obtained with the FIOOOW, F1280W, and 
F1800W filters. For each filter, the FASTR1 readout pattern was cho- 
sen with a two-point dither pattern. W1935 was observed with MIRI 
on 20 September 2022 using 15 groups per integration for FIOOOW, 
13 groups per integration for FIZ80W and 11 groups per integration 
for F1I800W. The total exposure time plus the overhead for the MIRI 
observing of W1935 was 1.03 h. W2220 was observed with MIRI on 18 
October 2022 using seven groups per integration for FlIOOOW, seven 
groups per integration for F1280W, and ten groups per integration 
for F1I800W. The total exposure time plus the overhead for the MIRI 
observation of W2220 was 0.54 h. 

As noted in ‘Sample’ above, both W2220 and W1935 have previously 
published photometry and W2220 has previously published spectra. 
Extended Data Table 1 lists all observables, both previous and new in 
this paper, included in the analysis that follows. 


Datareduction 

We used the official JWST science calibration pipeline (v.1.10.0) to 
reduce the NIRSpec G395H spectra starting from uncalibrated data 
downloaded from the Mikulski Archive for Space Telescopes (MAST). 
The pipeline comprises three separate stages. Stage 1 performs 
detector-level corrections (for example, bias subtraction, dark sub- 
traction and cosmic-ray detection) and ramp fitting to generate a 
count rate image for the individual uncalibrated image of each expo- 
sure. The resulting count rate images were calibrated by applying 
instrument-level and observing-mode corrections in Stage 2. Stage 
3 combines several calibrated exposures and extracts the spectrum. 
We optimized the aperture extraction location by using the relative slit 
position of the target to account for inaccuracies in the object coordi- 
nates and the celestial World Coordinate System. 

The flux uncertainties automatically propagated through the pipe- 
line were all null for the extracted spectrum due to the most recent 
reference flat files for NIRSpec having NaN values. We recalculated 
the flux errors for the extracted spectrum by combining in quadrature 
the Poisson variance (FLUX_VAR_POISSON) and read noise variance 
(FLUX_VAR_RNOISE) provided in the extracted spectrum file. 

Looking at the W1935 3.326 um spectral feature alone, we carefully 
examined each dither position to confirm that the feature was pre- 
sent. Although this is, in general, an area of the spectrum with a low 
S/N (average of 5-10 across the entire feature), we found the methane 
emission persisted in individual single exposure dithers, thus confirm- 
ing its presence. 

For both W2220 and W1935, we used the MAST-produced pipeline 
reductions of MIRI photometry, choosing the aper total vegamag 
column as our preferred magnitude. 


Radial velocity 

Given the high resolution of G395H data, we were able to compute radial 
velocities for both targets. All values reported are from a correlation 
with the models of ref. 44. Our values for W1935 and W2220, which are 
listed in Extended Data Table 1, are -36.9 + 5.1 and -53.2 + 2.8 km s”, 
respectively. We used these values in the banyan sigma kinematic code* 
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to check whether our targets belonged to any known moving group. 
For both objects, the full kinematics yielded a 99.9% field population 
probability. 

We also examined the total space velocity for each object and com- 
puted values of 42.02 + 5.33 and 55.33 + 2.82 kms‘ for W1935 and 
W2220, respectively. These values are consistent with the field brown 
dwarf population (see, for instance, ref. 46). Consequently, we have no 
evidence for youth or significant old age in the observables. Therefore, 
we choose a broad age range of 4.5 + 4.0 Gyr to estimate fundamental 
parameters in the SED approach (see below). 


SED construction and results 

Key to the analysis of this work was generating distance-calibrated 
SEDs inclusive of the newly obtained JWST data and all previously col- 
lected observables. To construct the SEDs, we used the open-source 
package SEDkit (ref. 14), which was first published and used by ref. 47 
to analyse the fundamental parameters of brown dwarfs. As described 
in ref. 47, we first combined the parallax with spectra and photometry. 
For W1935, the only spectrum available was the G395H data. For W2220, 
an HST grism spectrum was available from ref. 39 as well as the newly 
acquired JWST data. All photometry used in the SED for both objects 
is listed in Extended Data Table 1. 

Using SEDkit and the input data, we constructed the SEDs, integrated 
under the data as described in ref. 47 and calculated the L,,,, value for 
each object. To examine the similarities and differences for the two 
objects, Fig. 1shows the output SEDs for both W2220 and W1935 scaled 
to 10 pc. As described in the main text, the age estimate was paired 
with the evolutionary models of ref. 48 to acquire a radius range, and 
we then semi-empirically calculated 7,,,, log g and mass. All values are 
listed in Extended Data Table 1. 

Figure 1 shows the G395H portion of the SEDs overplotted. The SEDs 
are both excellent fits for each other but show important differences. 
The 3.326 um CH, bandis striking given that it is clearly inthe emission 
from W1935 but absent for W2220. There is also a difference in the 
3.8-4.3 um region, where W2220 appears to have the same molecular 
features (see Figs. 2 and 3 and Extended Data Figs. 4 and 5), but the gas 
is warmer. Hence, the flux is higher for W2220. which drives its bluer 
mid-infrared colour. 


Retrieval analysis 

We carried out retrieval analysis of the NIRSpec spectroscopy of our 
two targets using the Brewster package. Brewster is publicly available 
and was developed to model substellar atmospheres. It has been suc- 
cessfully applied to a range of brown dwarf atmospheres from cool T 
dwarfs to hot L-type subdwarfs, including cloudy and planetary mass 
objects 40, This is its first application in the Y dwarf regime. 


Retrieval method. Brewster consists of a forward model and sampler. 
The forward model produces a synthetic spectrum for a given set of 
atmospheric parameters. This is coupled to a Bayesian sampler, which 
explores the parameter space and estimates the posterior distribution 
for the forward model parameters given the data. 

The forward model is a one-dimensional, radiative-transfer, layered- 
atmosphere, model consisting of 64 layers distributed uniformly in 
log pressure space (in bars) between 2.3 and —4.0. The model evaluates 
the emergent flux using the two-stream source function technique of 
ref. 51, including scattering”. This requires each layer in the atmos- 
phere to be specified in terms of the wavelength-dependent optical 
depth, single scattering albedo and asymmetry parameter, as well as 
the temperature of the gas in the layer. The resultant spectrum is then 
processed to allow direct comparison to the data and the calculation 
of alikelihood for the fit. 

The parameters for the forward model fall into the following groups: 
gas-phase opacity, cloud opacity, temperature structure and global 
properties of the target. The gas-phase opacity is set by the choice of 


absorbing gases, their concentration and distribution in the model 
atmosphere. We included the following gases in this study: H,O, CH,, 
CO, CO,, NH;, H,S and PH,, In this case, we parameterized the concen- 
trations of the absorbing gases as vertically constant volume mixing 
ratios as in previous studies with this code. This approach provides 
substantial flexibility in arriving at possible solutions while still retain- 
ing computational simplicity. We neglected cloud opacity in this study 
and do not discuss this aspect of Brewster any further here but leave 
that work for a future study. 

We modelled the temperature structure using the same method 
as refs. 53,54. Briefly, we parameterized the temperature at 13 evenly 
spaced points in log pressure and evaluated the temperature in our 
atmospheric layers using spline interpolation from these. We avoided 
implausible discontinuities and so-called ringing by penalizing the 
second derivative of the thermal profile in the likelihood by using the 
following log prior on the temperature: 
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Inp(T)=- ay 


N 
È a- 27;+ T9- Zinn). (a) 
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This sums the discrete second derivative of the temperature Tat each 
level i, and weights it by y. The parameter y sets the degree to which 
the likelihood will be penalized by wiggles in the thermal profile and 
is included in our retrieval as a free parameter. By including this free 
parameter, the data can set the degree of smoothing imposed on the 
profile. We follow ref. 53 in adopting an inverse / distribution as the 
hyperprior for y, with properties specified in Extended Data Table 3. 

We parameterized the global properties of the target as follows. The 
radius of the source is encoded in the scale factor that is applied to the 
top-of-atmosphere flux produced by our radiative-transfer code to 
allow comparison to the flux received by NIRSpec. This scale factor is 
equal to R?/D?, where Ris the radius of the source and Dis its heliocentric 
distance. Our forward model uses R?/D’, and hence, this is the quantity 
that is estimated directly by the retrieval. We used our knowledge of the 
target’s distance to estimate the radius post hoc and incorporated the 
uncertainties in both the distance and the absolute flux calibration of 
the spectrum in the error estimate for the radius. The surface gravity 
for the source is parameterized as log g, where gis the gravitational 
acceleration at the altitude of our model atmosphere in cms”. We 
combined this parameter with our post hoc estimate for the radius 
to estimate the mass of the source. Finally, we included parameters 
to apply rotational broadening and radial velocity shifts to our model 
spectrum. We used the rotational broadening code provided by ref. 
55 to achieve the former. 

In addition to the free parameters described above, we hardcoded 
the following into our model atmosphere. In addition to the absorbing 
gases, we assumed that the atmosphere is dominated by H, gas, with 
He present at a solar H/He ratio. The total abundance of H,/He was set 
by the remaining fraction after the absorbing gas fractions have been 
accounted for. 

We calculated layer optical depths due to the absorbing gases using 
opacities sampled at a resolving power R = 30,000 taken from the com- 
pendium of refs. 56,57, with updated opacities described in ref. 15 using 
the same method as ref. 16. We generated the latest version of CH, 
line list” with broadening coefficients relevant to H,/He atmospheres 
using the computational methodology of ref. 59 and incorporated PH, 
opacities from refs. 48,60. We also included continuum opacities for 
H,-H, and H,-He collisionally induced absorption using the cross sec- 
tions from refs. 61,62 and Rayleigh scattering due to H,, He and CH, but 
neglected the remaining gases. We also included continuum opacities 
due to bound-free and free-free absorption by H (refs. 63,64) and 
free-free absorption by H; (refs. 65). 

The emergent spectrum was convolved with a Gaussian kernel of 
width 2.2 pixels to simulate the instrumental profile of the NIRSpec 
G395H. 


Inthis work, we coupled the forward model to the emcee sampler®, 
whichis the same method used in refs. 15-18,49. We used a log-likelihood 
function to assess the fit of the data to the model: 
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where the index irefers to the ith of n spectral flux points y, with errors 
S, which are compared to the forward model fluxes F; for the current 
parameter set x. To (at least partially) account for model and other 
unaccounted sources of uncertainties, we inflated our errors using a 
tolerance parameter, such that our data error s;is given by: 


s?=07+10°, (3) 


where o;is the measured error for the ith flux measurement and bis our 
tolerance parameter, which is retrieved® ©. The full set of parameters 
used and their priors are listed in Extended Data Table 3. 

We used 16 walkers per dimension. Following 10,000 iterations of 
burn-in, we ran our chains for blocks of 30,000 iterations, checking 
each time for convergence. In all cases convergence appeared to be 
achieved after the first 30,000 block, but we ran an additional 30,000 
in each case nonetheless. 


Retrieval results. We ran three models each for W1935 and W2220: 
(1) our baseline model as described earlier, (2) the baseline model 
adjusted with an additional prior such that T; < 7;,,, which excludes 
the possibility of a temperature inversion and (3) the baseline model 
adjusted to include H," as an absorbing gas at pressures less than 1 mbar. 
To distinguish the preferred model given our data, we calculated the 
Bayesian information criterion (BIC) for each case. For a set of mod- 
els, the one with the smallest (typically most negative) BIC will be the 
preferred model with the strength of the preference dependent on 
the value of the difference in the BICs, ABIC. By convention, the ‘win- 
ning’ model is defined as having ABIC = O and lower-ranked models 
have ABIC quoted relative to that. Reference 69 provided the following 
intervals for selecting between two models under the BIC: 

e 0< ABIC < 2: no preference worth mentioning; 

e 2< ABIC <6: positive; 

e 6 < ABIC < 10: strong; 

e 10 < ABIC: very strong. 


We found that for W1935, the baseline model is very strongly pre- 
ferred over the model that includes H;*. Both of these models are very 
strongly preferred over the model that rules out a temperature inver- 
sion. By contrast, for W2220 we found that there is no preference worth 
mentioning between models with and without a temperature inversion, 
which are both strongly preferred over the model that includes H,”. This 
difference arises because W1935 shows emission in the CH, q-branch 
at 3.326 um whereas W2220 does not. The ABIC values for each model 
are given in Extended Data Table 4. 

Note that the inclusion of H,” in the model had no impact on the 
quality of the fit, the maximum likelihood or the values of the other 
parameters. The abundance of H,” was an upper limit oflog fy =-6. 
We, thus, conclude that H;* is not detected in either source. 

Parameter estimates (not related to the temperature structure) for 
W2220 are indistinguishable between the models that allow for a tem- 
perature inversion and those that do not. A comparison of the posterior 
distributions for models with and without a temperature inversion is 
shown in a corner plot in Extended Data Fig. 2. 

The thermal profiles are also extremely similar, with only nonsigni- 
ficant differences between the two. The two retrieved thermal profiles 
are compared in Extended Data Fig. 3. The two profiles are identical at 
pressures deeper than log[P (bar)] =- 2.0 and do not deviate signifi- 
cantly from one another at shallower pressures. Extended Data Fig. 4 


shows that the retrieved model spectra are also indistinguishable 
anda similarly good fit to the data as the baseline model for W1935 
(Fig. 2). 

The posteriors for the non-temperature-related parameters for 
W1935 are also very similar between the two models, as shown in 
Extended Data Fig. 5. Although there are no significant differences 
between the compositions, there is a marginal trend across all the 
absorbing gases towards higher abundances in the baseline model 
(which allows for a temperature inversion). 

Figures 2 and 3 demonstrate clearly that the model that allows a 
temperature inversion is the only one that is able to fit the CH, emis- 
sion feature at 3.326 um, which is absent in the spectrum of W2220. 

Our retrieval results for both objects show some differences from 
those inferred from our SED-based luminosity and evolutionary model 
predictions (Extended Data Figs. 5 and 2 and Extended Data Table 1). 
In both cases, the masses implied by our retrievals are slightly higher 
than the lo upper limit suggested by evolutionary models. The dif- 
ference for W1935 is 1.30, and for W2220, it is 1.10. It is not unusual for 
retrieval analyses to disagree with evolutionary models’ predictions 
of logg, radius and mass, particularly in the absence of any strong prior 
evidence to provide empirical constraints®"*°°”°, In addition, these 
retrievals cover a relatively narrow wavelength range that lacks recog- 
nized gravity-sensitive spectral features. Establishing the presence and 
nature of any possible biases to higher gravity and mass estimates is 
beyond the scope of this work and does not impact our central results 
or conclusions. 

Extended Data Fig. 6 compares the contribution functions for maxi- 
mum likelihood retrieval models for W1935 and W2220. The contribu- 
tion function in an atmospheric layer lying between pressures P, and 
P, is defined as: 


BA, T(P) J? dt 


C(A, P) = (4) 


C(A, P) effectively maps the depth in the atmosphere from which the 
flux observed at a given wavelength originates. Extended Data Fig. 6 
demonstrates that the CH, emission seen in W1935 originates from 
pressures shallower than 0.01 bar. Our model finds gas some 300 K 
hotter than is retrieved when no inversion is permitted. 

For simplicity, we used only cloud-free retrieval models. It has 
been well documented that the omission of clouds in a retrieval of 
an atmosphere that contains clouds can bias the thermal profile to a 
more isothermal gradient, as the retrieval mimics the effect of clouds 
in screening hotter deeper layers in the opacity windows between 
molecular absorption features»””°. Hence, the kinks seen in the 
retrieved thermal profiles may be suggestive that clouds should be 
considered for future retrieval studies in this temperature regime. 


Data availability 

TheJWST data in this paper are part of GO programme 2124 (principal 
investigator:J.K.F.) and are publicly available from MAST (archive.stsci. 
edu/) under that programme ID. The HST WFC3 spectrum of W2220 
is available from ui.adsabs.harvard.edu/abs/2015Ap]J...804...92S/ 
abstract. 


Code availability 

The data reduction pipeline jwst can be found at jwst-pipeline.readthe- 
docs.io/en/latest/. The Brewster code is open source and available at 
the following GitHub repository: github.com/substellar/brewster. 
Similarly, the SEDkit code is open source and available at github.com/ 
hover2pi/sedkit. The set-up that yields the results presented herein is 
discussed in the Methods. 
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Extended Data Fig. 1| The Spitzer color-magnitude diagram for cold brown unceratinty on both color and absolute magnitude is plotted. We also label and 
dwarfs. The sources in program 2124 are green filled squares. Both W2220 and highlight as a black filled circle the coldest known brown dwarf 0855-0714 to 
W1935 are highlighted by a five-point star to emphasize their position. The lo demonstrate the extreme of Y dwarf parameters. 
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Extended Data Fig. 2| The retrieved composition of W2220, along with gravity, and derived quantities for mass, radius, metallicity and C/O. Posteriors for 
models with and without a temperature inversion are shown in purple and dark cyan respectively. 
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Extended Data Fig. 3 | The retrieved thermal profiles for W2220 with and 
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Extended Data Fig. 4| The JWST G395H spectrum for W2220. Overlaidin The 67% confidence intervals in the model posteriors are shaded in dark cyan 
purple and dark cyan are the median retrieval models for the source with and purple, but are of comparable extent to the line width. 
and without temperature inversions. Data uncertainties are shaded in grey. 
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Extended Data Fig. 5| The retrieved composition of W1935, along with gravity, and derived quantities for mass, radius, metallicity and C/O. Posteriors for 
models with and without a temperature inversion are shown in purple and dark cyan respectively. 
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Extended Data Fig. 6 | The contribution functions for W1935 and W2220 maximum likelihood retrieval models. The CH, emission seen in W1935 originates 


from pressures shallower than 0.01 bar. 


Extended Data Table 1| Parameters of interest for W1935 and W2220 


Astrometry 


a (J2000) 
6 (J2000) 
W 
Distance 
Ha 

Hô 


Kinematics 


X 
Y 
Z 
U 
V 
W 
RV 


Total Velocity 
Photometry 


WISE W1 
WISE W2 
Spitzer [3.6] 
Spitzer [4.5] 
JWST F1000W 
JWST F1280W 
JWST F1800W 


Fundamental Parameters 
Age 
Teff 


log(Lpoi/Lo) 
Radius 
Mass 


log g 


Data from refs. 7,36,37,39,71. 


W1935 


293.827684 
—15.772363 
69.343.8 
14.43+0.79 
290.2+11.6 
43.1411.5 


12.69 +0.69 
5.49+ 0.30 
-4.14 + 0.22 
-41.20+ 4.52 
-6.03+2.12 
-5.61+ 1.87 
-36.9+5.1 
42.025.33 


>Y1 

23.930.33 
18.75+0.30 
15.79+0.06 


18.5140.03 
15.530.02 


13.740+0.005 
13.126+0.007 
12.107+0.017 


4.5+4.0 
482+38 
-6.30.1 
0.950.14 
6-35 
4.70.5 


W2220 


335.231035 
—36.471677 
95.5+2.1 
10.47+0.23 
290.1+0.9 
-97.1+0.9 


5.67 + 0.12 
0.70 + 0.02 
-8.77 + 0.19 
-39.6141.54 
-11.10+0.25 
37.00+2.35 
-53.2+2.8 
55.3342.82 


YO 

20.91+0.09 
20.640.05 
20.96+0.08 
21.330.15 
18.29+0.15 
14.81+0.02 
17.20+0.06 
14.74+0.02 


13.354+0.006 
12.979+0.009 
12.276+0.021 


4.5+4.0 
480+41 
-6.4+0.1 
0.94+0.14 
6 — 35 
4.7+0.5 


ref 


[7], [7] 
[7], [7] 
[7], [7] 
[7], [7] 
[7], [7] 
[7], [7] 


This paper 
This paper 
This paper 
This paper 
This paper 
This paper 
This paper 
This paper 


[7], [39] 
[71] 

[37], [71] 
[71] 

[71] 

[7], [7] 

[7], [7] 
[36], [7] 
[36], [7] 
This study 
This study 
This study 


This study 
This study 
This study 
This study 
This study 
This study 
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Extended Data Table 2 | Retrieved Gas Abundances and Fundamental Parameters 


log fgas 
H20 
CHa 
CO2 
CO 
NH3 
HeS 
PH3 


Metallicity and C/O 


[M/H] 
C/O 


W1935 


+0.07 
za tet 
rostit 
—7.65 0.11 


+0.08 
—5.09 011 


+0:08 
—3.97 oii 


W2220 
+0.10 
>o 
—4.92 017 


Extended Data Table 3 | Parameters and priors adopted for the retrieval analysis 


Parameter Prior 

gas fraction, (Xgas) log-uniform, log Xgas > —12.0, De as Xgas < 1.0 

thermal profile, T uniform, constrained by 0.0 K < T; < 5000.0 K 

profile smoothing parameter, y uniform, 0 < y < 10° 
gravity, log g uniform, constrained by 1Mjup < gR?/G < 80M Jup 
scale factor, R?/D? uniform, constrained by 0.5Rjup < R < 2.0Rjup 
rotational velocity, v sin t uniform, 0 kms~! < vsini < 100 kms7! 

radial velocity, Vrad uniform, —250 kms~! < vraq < 250 kms~! 


tolerance factor, b uniform, log(0.01 x min(o?)) < b < log(100 x maz(a?)) 
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Extended Data Table 4 | ABIC values for the models tested for each target 


Model ABIC 
W1935 

Baseline (allowed inversion) 0 
No inversion 63.6 
Baseline With Hy 9.9 
W2220 

Baseline (allowed inversion) 1.8 
No inversion 0.0 


Baseline With Hj 9.9 


